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Abstract 

Ultracold polar molecules trapped on an optical lattice is a many-body system 
that, under appropriate conditions, may support collective excitations reminiscent 
of excitons in solid state crystals. Here, we discuss the rotational excitations of 
molecules on an optical lattice leading to rotational Frenkel excitons. Apart from 
solid hydrogen, there is no other natural system that exhibits rotational excitons. The 
rotational excitons have unique properties that can be exploited for tuning non-linear 
exciton interactions and exciton-impurity scattering by applying an external electric 
field. We show that this can be used to explore the competing role of the dynamical 
and kinematic exciton-exciton interactions in excitonic energy transfer and to study 
Q-f quantum localization in a dynamically tunable disordered potential. The rotational 

excitons can also be used as a basis for quantum simulation of condensed matter 
models that cannot be realized with ultracold atoms. As an example, we discuss 
the possibility of engineering the Holstein model with polar molecules on an optical 
lattice. 

1 Introduction 

Ultracold atoms offer the possibility to study few- and many-body quantum systems with 
exquisite control over microscopic interactions. This has led to spectacular experiments pQ , 
bridging different areas of physics. While the research directions stimulated by experiments 
with ultracold atoms are very diverse, one can identify several overarching trends that 
gained momentum in recent years. The most prominent of these trends is the effort aimed at 
understanding quantum phase transitions [2] , leading to the experimental studies of Bose- 
Einstein condensation (BEC) [3], bosonic superfluidity pi], quantum magnetism [5], many- 
body spin dynamics [6], Efimov states [TJ, Bardeen-Cooper-Schrieffer (BCS) superfluidity 
[8] and the BEC-BCS crossover [9]. With the development of techniques for trapping 
atoms in periodic potentials of optical lattices [10] and single atom detection [11] . ultracold 
atoms became an ideal platform for quantum simulation of lattice models used in solid 
state physics p2j. As demonstrated by the realization of the Mott insulator - superfluid 
transition with atoms trapped in a three-dimensional optical lattice [13] , the experiments 
with ultracold atoms hold the promise of insight into the details of the Hubbard model, 
which may help unravel the mechanism of high-T c superconductivity 



1 



The development of experimental methods for the production of ultracold polar molecules 
[T5] has widened the possibilities for quantum simulation of condensed matter models to a 
great extent. A combination of the rotational, spin and hyperfme degrees of freedom with 
the long-range dipolar interactions enabled by the dipole moment of molecules allows for 
engineering a great variety of lattice models that cannot be realized with atoms [THl fT7] . 
Of particular interest is the possibility of creating quantum phases with topological order 
|18j . which are resilient to perturbations preserving the topology and are, therefore, ideal 
for quantum computation. To this end, a major effort of current experiments is focused on 
preparing a dense ensemble of ultracold polar molecules trapped in optical lattices [19]. 

Another focus of research with quantum degenerate gases is on emergent phenomena, 
such as solitons [20] . rotons [21], vortices [22] , spin waves [23] and polarons [21]. These 
experiments aim to elucidate emergence in natural systems, and may potentially lead 
to the development of novel ultra-sensitive sensors of gravity and electromagnetic fields. 
While these studies cover a wide range of collective dynamics of ultracold atoms and 
molecules in a single quantum state, much less is known about the effect of internal degrees 
of freedom of ultracold particles. For example, the role of rotational transitions in the 
excitation spectrum of a molecular Bose-Einstein condensate had not been addressed until 
very recently [25]. Yet, the internal degrees of freedom can be used to explore new regimes 
of collective phenomena, especially in molecular systems that provide a dense spectrum of 
internal excitations [261 127] . 

In this article, we consider rotational excitations of polar molecules trapped on an 
optical lattice. These excitations give rise to rotational excitons analogous to collective 
electronic excitations in molecular crystals. In contrast to excitons in natural solids, the 
properties of rotational excitons can be dynamically controlled by tuning the energy level 
structure of the trapped molecules, which can be used to study new regimes of Frenkel 
exciton physics not accessible in natural solid state crystals. For example, we show that 
non-linear interactions of excitons can be tuned to examine the competition between the 
dynamical and kinematic effects and demonstrate that rotational excitons can be used to 
study quantum localization in a dynamically tunable disordered potential. The rotational 
excitons can also be used as a basis for quantum simulation of condensed-matter models 
that cannot be realized with ultracold atoms. As an example, we discuss the possibility 
of engineering the Holstein model with polar molecules on an optical lattice. In order to 
present these results in the context of current research, we briefly describe the related work 
on quantum simulation of many-body Hamiltonians with ultracold atoms and molecules 
and the basics of Frenkel exciton physics. 

2 Quantum simulation of lattice models 

The difficulty of simulating a quantum many-body system on a classical computer increases 
exponentially with the number of quantum states. An alternative, currently at the focus of 
detailed research, is quantum simulation [28]. Quantum simulation involves the design of 
a controllable quantum system in order to simulate the properties of another, more com- 
plicated system (analog quantum simulator). This idea is generally attributed to Feynman 
[29] . although quantum simulation is already mentioned in an earlier publication by Manin 
[30] . The key ingredient of quantum simulation is the mapping of the Hamiltonian of 
the simulated system onto the Hamiltonian of the simulator with controllable parameters. 
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Tuning the parameters of the controlled system can then be used to map out the phase 
diagram of the simulated system. The recent literature offers many ideas about how to 
build quantum simulators based on ultracold quantum gases [31], ultracold trapped ions 
[32] . single photon sources and detectors [33], ultracold Rydberg atoms [31], circuit QED 
[35] , quantum dots and N-V centers in diamond [36J. 

Ultracold atoms trapped in a periodic potential of overlapped laser beams (optical 
lattices) is considered as one of the most promising systems for quantum simulation of 
lattice models [101 |37]. The geometry and strength of the optical lattice potentials can 
be controlled by adjusting the laser intensity and the beam overlap angle, which enables 
control over the translational motion of atoms in the lattice [38J. The separation of the 
lattice sites in an optical lattice is equal to half the wavelength of the trapping laser 
field. Recent technological breakthroughs permit the detection of ultracold atoms with 
sub- wavelength resolution, allowing for single-site addressing [39] . 

2.1 Lattice models 

The most prominent example of quantum simulation with ultracold atoms is the realization 
of the Hubbard model gQJ [T3J 

H H = ~t + -n H ™*a™".t + £ E*n n ,a, (1) 

where a =t or J,, & n a and c n cr are the creation and annihilation operators for fermions, 
satisfying the anti-commutation rule {& n CT , c m>a } = 6 nm , n„ )(7 = & n a c n ^ is the particle 
number operator, and E a is the energy of particles in state a. The angle brackets (n, m) 
indicate that the summation is over nearest neighbors only. In this simple form, the 
Hubbard Hamiltonian accounts for nearest neighbor tunnelling (t) and on-site interactions 
for fermions in different spin states (U). The fermions with the same spin experience hard- 
core repulsion. This model is widely used for the studies of high-T c superconductivity [JT] 
and quantum magnetism [4"2"] . 

With bosonic atoms, it is possible to realize the Bose-Hubbard model: 

H B h = -* Ofim + — h n (h n - 1) + E Y n n , (2) 

(n,m) ti n 

where b\ and b n are the creation and annihilation operators for bosons, satisfying the com- 
mutation rule b m ] = 8 nm , and h n = b}j) n is the particle number operator. This model 
has been used in many studies of the properties of bosonic gases such as superfluidity [43] . 
If U <f, the particles are delocalized around the lattice and the ground state of the sys- 
tem is a superfluid. This Hamiltonian also allows for modelling the properties of fermionic 
systems when quantum statistics plays no role and the properties of the system are entirely 
determined by the relative efficiency of hopping (governed by t) and the interaction (U). 
In particular, when U 3> t (the on-site repulsion dominates over the inter-site hopping), 
this bosonic Hamiltonian reproduces transition to the Mott insulator phase [Sj, where 
the number of atoms per lattice site is fixed. This transition is typical for electrons in a 
metal. It is thus possible to realize the same physical phenomenon with particles satisfying 
different quantum statistics. We shall discuss this in detail in section 12.21 
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Recent experiments demonstrated the possibility of creating the Mott insulator phase 
with ultracold atoms filling up to 99% of the lattice sites |45j. A major thrust of current 
research is to extend these experiments to ultracold molecules. The experiments with polar 
molecules in optical lattices offer the possibility of realizing lattice models with long-range 
interactions [TCI H7J H5J H7]. In the Bose-Hubbard model the inter-site interactions are 
accounted for by adding to the Hamiltonian (J2]) a term H v = K^m^n^m describing 

(n,m) 

density-density correlations. In the limit of the strong on-site repulsion, when (n n ) = or 
1, the system reduces to the t-V model, described by the Hamiltonian 

H t -v = -t Hfim + E J2n n + V ^ n n h m (3) 

(n,m) n (n,m) 

with the constraint 

St£t|$) = IXWi = o (4) 

for any physical state of the crystal |$), which accounts for the infinite on-site repulsion. 
With V = this Hamiltonian corresponds to the lattice analog of the To nks- Girardeau gas 
[HISS]- This model is equivalent to an anisotropic spin- 1/2 XXZ model [50] . 

In addition, molecules allow for the possibility of realizing a variety of lattice spin 
models [T2], such as 

H* = ~J2 ( J *« + J y« + Jz«) <> ( 5 ) 

(n,m) n 

where a^ ,v,z are the Pauli matrices for lattice site n. This is the so-called Heisenberg 
model. If the interactions are isotropic in the plane perpendicular to the z axis, i.e. for 
Jx — Jy, this Hamiltonian reduces to the XXZ^model [5T]. In the isotropic Heisenberg 
model, J x = J y = J z . The Xy-model is for the case when J z = 0. The scalar version of the 
Heisenberg model ( J x = J y = 0) is the Ising model. In the Ising model, spins are treated 
as scalars, which can take one of two values: s n = ±1. The Hamiltonian becomes 

Hi = - JnmSnSm ~ S n . (6) 

n,m n 

The Ising model can be realized with ultracold molecules in the spin-less X S electronic 
state trapped on an optical lattice in a Mott insulator phase with one molecule per lattice 
site. The rotational levels of X S molecules in an external electric field form an isolated 
two- level system, illustrated in Figure [Th. The ground \g) and excited |e) states can be 
used as the spin states s n . The coupling constant J nm is determined by the dipole-dipole 
interaction between molecules in different lattice sites. The rotational excitation \g) — > \e) 
thus leads to spin waves of the Ising model. In order to simulate more complicated spin 
models, such as Eq. (jSJ, it is necessary to use molecules with more complex structure 

2.2 Bosons, Fermions and Paulions 

The Hamiltonians presented in section 12.11 describe very different physical systems, con- 
sisting of particles with all possible statistics (bosons, fermions, spin matrices and pseu- 
dospins). They are, however, mutually related. Some of these Hamiltonians can be mapped 
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onto each other, even when the quantum statistics of mapping and mapped particles are 
different. 

Consider an ensemble of arbitrary two-level systems on a lattice. The state of a two- 
level system in lattice site n can be characterized by the operators p\ and p n describing, 
respectively, the creation and destruction of the excited state. Operators for different 
lattice sites necessarily commute, as they act on different variables. On the other hand, 
the same site cannot accommodate more than one excitation, which is a feature of Fermi 
statistics. The commutation relations for the operators p n are 

PnPln ~ PlnPn = (u ^ m) 

(7) 

PnPn + PlPn = 1- 

The operators with such "mixed" statistics are called paulions (in condensed matter 
physics [52]) or hard-core bosons (in atomic physics [IH]). Using the p-operators, we can 
write the following general Hamiltonian 

A / / 

,mP\iP\nPnPrn (8) 

n n,m n,m 

which is identical to the hard-core Hamiltonian ([3]) with the constraint (J4]) absorbed into 
the statistical properties of paulions ([7]). The prime over the sum symbols indicates that 
n m. For V = 0, i.e. when the density-density correlations are absent, the Hamiltonian 
OH]) is identical to the Ising Hamiltonian fl6]) with h = 0. This Hamiltonian is also used to 
model Frenkel excitons in solid state molecular crystals [53J. 

Any unitary transformation preserves the commutation properties of the bosonic and 
fermionic operators. For example, the Fourier transform of bosons must produce bosons. 
Therefore, many articles have been devoted to mapping paulions onto "effective" particles 
with bosonic or fermionic statistics^. Girardeau showed [H] that in ID the many-body 
wave function of hard-core bosons - paulions - corresponds (up to a sign) to a many-body 
wave function of a gas of fictitious non-interacting spin-less fermions. The exact map- 
ping was proposed later by Chestnut and Suna for a ID system in the nearest neighbor 
approximation [51]. This is a variant of the Jordan- Wigner transformation [51], and it 
works very well in ID. However, in higher dimensions it may cause problems: the effective 
fermionic operators are non-local, i.e. the expression for the fermionic operators in site 
n is dependent on the occupation numbers at other sites. That is why this transforma- 
tion is not as effective in dimensions higher than one. Another possibility is to use the 
Agranovich-Toshich transformation [57], which expresses p-operators through an infinite 
series of bosonic operators b n and 6^, in a way ensuring that for any number of bosons per 
lattice site N^> = 6^6 ra = 0,1,2,3... the eigenvalues of the number operator for paulions 
N n = p\j)n are only 1 or 0, so that the unphysical states with N n > 1 do not occur. 
Then the paulionic Hamiltonian (jSJ) reduces to a sum of terms consisting of the same 
number of creation and annihilation bosonic operators, which we schematically denote by 
H v = a ubl ll ---bl lu b mi ...b mi/ ; v = l...oo. Keeping the pairwise interactions only and 

{n,m} 

1 We note that similar transformations exist also for spin matrices er+ = + ia~{, a~ — a® — ia^ and 
a z n . They can be mapped onto fermionic operators by the Jordan- Wigner transformation [54] . and onto 
bosonic operators by the Holstein-Primakoff transformation |55j . 
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assuming that \t\ -C E , one obtains: 

H = E Y, P n k + E tn,nMfm + Y! V rhmKHJ>„b m + # kin , (9) 
n n n,m 

where 

N 

H kiQ = -2E Y blPXL (10) 

n=l 

Thus, the Agranovich-Toshich transformation reduces the pauilonic (hard-core boson) 
Hamiltonian to a Hamiltonian describing a gas of bosons interacting via delta-like pairwise 
interaction with the strength 2Eq. Higher-order corrections to the interaction energy can 
be obtained by including the omitted terms. The Hamiltonian describes the kinematic 
interaction, which we discuss in section 14.11 In contrast to the Wigner- Jordan transforma- 
tion, this approach is very effective in 3D, where can be treated as a perturbation. For 
2D systems, the delta-like scattering with the magnitude 2Eq, which is the largest energy 
scale of the problem, is very strong and does not allow for a perturbative treatment [58] . 
In ID, the effect of scattering is even stronger and the transformation to bosons makes no 
sense at all. To summarize, in ID paulions are well described by a gas of fermions, in 3D 
by bosons, in 2D they are something in between. Finally, we note that in principle it is not 
necessary to use any of these transformations; one can work directly with paulions taking 
into account their commutation relations (J7J). 

3 Rotational Frenkel excitons 

While the majority of research with atoms and molecules on optical lattices has so far 
focused on quantum transport of particles in a lattice potential or quantum simulation 
of lattice models, such as the ones described above, the ability to control the structure 
of molecules trapped on an optical lattice can also be exploited to study quantum en- 
ergy transport and collective excitations reminiscent of excitons in solid state crystals. In 
this section, we discuss the formation of rotational excitons. These excitons have unique 
properties and can, in turn, be used as a basis for quantum simulation of new physical phe- 
nomena and new lattice models that cannot be realized with atoms on an optical lattice. 
This is discussed in sections 0] - |6j 

Rotational excitation of polar molecules trapped on an optical lattice in the Mott in- 
sulator phase with one molecule per site gives rise to the formation of rotational Frenkel 
excitons [59]. The Frenkel exciton is a charge-less quasiparticle, which describes the ex- 
citation transfer in molecular crystals [60]. Here, we consider the transition between the 
absolute ground state \g) of trapped polar molecules and their first rotational excited state. 
We assume that the molecules reside in the ground vibrational state of the ground elec- 
tronic state X S. In the presence of a dc electric field £/, the three-fold degeneracy of the 
rotationally excited state \N = 1, Mjv) is lifted and the state with the projection = 
of the total angular momentum is detuned from the states with = ±1 (Figure [T^i). 
We assume that the detuning is large enough so that the latter can be disregarded and 
consider the isolated two-level system of the state \g) and the rotational excited state with 
Mat = 0, denoted by |e). 
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The molecular states are the eigenstates of the Hamiltonian 

Ht° l) = B e N 2 n - d n ■ S f , (11) 

where B e is the rotational constant, N n is the operator of the angular momentum, and 
d n is the dipole moment of the molecule in site n. The field-dressed states are linear 
combinations of the field-free rotational states 

|ff) = Ea#/) \N,M N = 0), 

N 

(12) 

\e)=EM£f) \N,M N = 0). 

N 

The coefficients and 0n are determined by the electric field strength Sf. In the limit 
Sf —> 0, the states \g) and |e) become, respectively, \N = 0, M N = 0) and \N = 1, M N = 0). 

It is convenient to introduce the transition operators P n defined by the equations 
Pn\9m) = S nm \e n ) and P n \e m ) = S nm \g n ). As discussed in section 12.21 these operators 
describe hard-core bosons, or paulions, and the excitonic Hamiltonian describes a hard- 
core boson gas with long-range interactions, which is equivalent to the Ising model. The 
rotational excitons can thus be mapped onto other systems described by the t—V model, 
or the Ising model. 

Using the P n -operators, the total Hamiltonian for M molecules on an optical lattice 
can be written as [53] 

Htot = Hexc + H dyn + H n _ c , (13) 
where the first term accounts for the excitation transfer between lattice sites 

^cxc = E EoHPn + Y,' J ( n - ™)PlPm (14) 

n n,m 

with the constant J(n — m) determined by the matrix elements of the inter-molecular 
excitation transfer due to dipole-dipole interaction 

J(n - m) = (e n g m \V dd (n - m)\g n e m ). (15) 

The second term describes non-linear interactions between excitons: 

^dyn = D ( n ~ m)PlPiP n P m , (16) 

n,m 

with the interaction constant D(n — m) determined by 

D(n -m) = (e n e m \V dd (n - m)\e n e m ) + 

(g n gm\V dd (n - m)\g n g m ) - 2(e n g m \V dd (n - m)\e n g m ). (17) 

Finally, H n _ c contains terms that do not conserve the number of molecular excitations in 
the system: 

H n _ c = (e n e m \V dd (n - m)\e n g m ) (P n + Pty + 

n,m 

1 ' * 

^E ( e ne m \V dd (n - m)\g n g m )[P n P m + P^Pj n ) + 

n,m 

/ r 

+E [(9n9m\V dd (n - m)\g n e m ) - 

n,m 

(e n e m \V dd {n - m)\e n g m )] (P ] n P n P m + P^P^Pn)- (18) 
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If \g) and |e) are states of well-defined parity (such as, for example, the rotational 
states at zero electric field), the matrix elements D(n — m) as well as the linear and cubic 
terms in Eq. (fl~8l) must vanish. This is the case for molecular solids, such as anthracene 
and naphtalene, which are often considered as prototype systems for (electronic) Frenkel 
excitonsd. For an ensemble of molecules in an optical lattice, the magnitude of these terms 
can be tuned by applying an external electric field, which breaks the inversion symmetry 
(parity) of the molecules. To our knowledge, the role of the odd-P terms in Eq. f|T8|) has 
not been explicitly studied and remains an interesting open problem. The possibility to 
tune these terms by varying the field dressing of the states \g) and |e) may lead to new 
interesting phenomena. In what follows, we assume that the effect of H n _ c is small and 
neglect this term. 

Due to the translational invariance of the optical lattice, the Hamiltonian H exc can be 
diagonalized by the transformation: 

^ = ^Ee-^Pt(g), 

(19) 

to yield 

H cxc = Y, E ( ( i) p Kq)P(q), (20) 

where P^(q) and P{q) describe the Bloch plane waves, q is the (linear) momentum, and 
E(q) = Eq + J2 J(n)e iqn is the exciton energy forming a quasi-continuous band. For M 3> 1, 

n 

we can consider the wave vector q and the energy E(q) to be continuous variables. 
In the nearest neighbor approximation, 

E NNA (q) = E + 2J cosaq (21) 

and 

J = (e n gn+i\Vdd(a)\gnen+i) , (22) 

where a is the lattice constant. The energies are shown in Figure [T] for different angles 
9 between the one-dimensional molecular ensemble and the electric field £f (which deter- 
mines the sign of J and, consequently, the shape of the dispersion curve). In the nearest 
neighbor approximation, A = 4 J is the exciton bandwidth. The magnitude of J is a quan- 
titative measure of the collective excitation effects (in particular, h/J is the timescale of 
the excitation transfer between molecules in adjacent lattice sites). Note that J = at 
9 = 9* = arccos(l/v^3): at this angle the dipole-dipole interaction vanishes, and (in the 
dipole approximation) the molecules become decoupled. When this happens, no excitation 
transfer dynamics can occur. For polar molecules with permanent dipole moments of a few 

2 The terms of the kind PP and P'P', which do not conserve the number of particles, do not vanish in 
natural molecular crystals with the inversion symmetry. The usual model neglecting these terms is known as 
the Heitler-London approximation |53j . These terms can be accounted for as corrections, whose magnitude 
depends on the ratio between the dipole-dipole interaction matrix element J (a) and the excitation energy 
of the molecules Eq, as well as on the structure of the molecules [61] . 
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Debye, the values of J can be tens of kHz (for LiCs on an optical lattice with a = 400 nm, 
J ~ 25 kHz). 

For small q, where the dispersion is parabolic, we can introduce the effective mass m e g = 
—h 2 /Ja 2 . For 9 > 9* (9 < 9*) the constant J is positive (negative), and, consequently, 
the effective mass is negative (positive). Quasiparticles with negative effective mass have 
counterintuitive propagation properties: they move in the direction opposite to their wave 
vector (negative refraction [62]). Thus, by changing the angle 9 it may be possible to tune 
the excitons from the regime of normal propagation to the regime of negative refraction. 




m P ff< 





Figure 1: Panel (a): Rotational energy levels of a X E molecule in the presence of an electric 
field. The dashed lines show the positions of the field-free states. Panels (b) and (c): 
Energy of rotational excitons in a one- dimensional array of LiCs molecules trapped on an 
optical lattice with the lattice constant a = 400 nm in the presence of an electric field 
of 1 kV/cm. Full curves - exact numerical calculation; dashed curves - analytical result 
of Eq. (f2~Tl) . The electric field is directed perpendicular to the array axis (panel b) and 
parallel to the array axis (panel c). 

The possibility of tuning the parameters of the exciton Hamiltonians ( fl4"l) , ( fl6l) and ( |T8l) 
in an ensemble of polar molecules on an optical lattice opens the possibility to access new 
regimes of Frenkel exciton physics that cannot be observed in solid state crystals. Here, 
we discuss a few examples that we presented in our recent papers [59], [63], [64"] . 



4 Nonlinear interactions of Frenkel excitons 

In this section we study non-linear interactions between rotational Frenkel excitons and 
show that they can be dynamically tuned by an external electric field. This is important 
for applications of excitons in quantum information processing [66J, where excitons can 
be used as qubits. In particular, we consider the interplay of the dynamical interactions 
arising from Eq. f fl6|) and the kinematic interaction arising from the hard-core boson nature 
of the exciton operators and show that the relative importance of these two interactions 
can be tuned. 

4.1 Dynamical and kinematic interactions 

There are two types of non-linear interactions for Frenkel excitons. The first is the dynam- 
ical interaction given by Eq. (fl6|) . This interaction is determined by the matrix elements 
D[n) of the dipole-dipole interaction (see Eq. (fiT|) ). In the wave vector representation, the 
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dynamical interaction describes wave vector conserving scattering between two excitons 
with the exchange of momentum q 

^ = Jf E D(q)P\k 1 + q)P\k 2 -q)P(k 1 )P(k 2 ). (23) 

ki,k 2 ,q 

where D(q) = J2 D{n)e~ iqn . 

n=\ 

The second interaction mechanism originates from the hard-core repulsion of molecular 
excitations. Since excitons are paulions, or hard-core bosons with long-range interactions, 
they effectively repel each other when placed in the same lattice site. As a result, the direct 
products of P(A;)-operators are not the eigenstates of the exciton Hamiltonian H exc (fT4l . 
as they would be for bosons or fermions. The commutation relations for paulions ([7]) yield 
for a two-exciton state |$(&i, ^2)) = \P (ki)P (^2)) the following equation: 

# exc mh,k 2 )) = [E(h) + E{k 2 )\ fc 2 )>+ 

1 (24) 
+- Y, [E(qi)+E(q 2 )]\<S>( qi ,q 2 )). 

J* 91.92 

91+92 = fc l+ fc 2 

The second term on the right-hand side may be interpreted as total wave vector conserving 
scattering between plane- wave-like one-exciton states |-P^(<?i)) [52J. This type of scattering 
results from the kinematic interaction [67]. This interaction can also be included as a 
perturbative term H^ lTL in the total Hamiltonian for the effective bosons introduced via the 
Agranovich-Toshich transformation (section 12.21) . The non-linear properties of excitons 
(in the two-body interaction approximation) are determined by the balance between the 
dynamical and kinematic interactions. 



4.2 External field control of exciton— exciton interactions 

In the presence of an external electric field and with account of the dynamical interaction, 
the rotational Frenkel excitons are described by the Hamiltonian H n \ = H exc + Hd y n given in 
Eqs. (THj) and ffT6j) . In the nearest neighbor approximation H exc and H^ yQ are parametrized 
by two constants: D = D(a) and J = J (a). Since J determines the propagation prop- 
erties of excitons, it is implicitly related to the strength of the exciton-exciton kinematic 
interaction: roughly speaking, the strength of repulsive interaction (scattering) between 
two excitons at the same lattice site is proportional to \J\. The constant D describes the 
dynamical interaction. The paulionic corrections do not affect Hdyn so one can use the 
bosonic commutation relations for the P-operators in this term, both in the site and wave 
vector representations. 

Both D and J can be controlled by varying the magnitude and orientation of the electric 
field (see Figure EJ). Panel (a) of Figure [2] shows the dependence of D and J for a one- 
dimensional array of LiCs molecules on the field magnitude at a fixed angle between the 
array axis and the direction of the electric field. Note that D vanishes at zero electric field 
(see Eq. (fTTt) ). For £f = 0, the field-dressed states \g) and |e) reduce adiabatically to the 
pure rotational states: \g) — > \N = 0, = 0) and |e) — > \N = l,M^ = 0), and the 
matrix elements determining D vanish according to the selection rules. Panel (b) shows 
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Figure 2: (a) D and J as functions of electric field magnitude at a fixed angle 9 = 90°. (b) 
D and J as functions of 9 at a fixed electric field Sf = 6 kV/cm. The calculations are for 
a ID ensemble of LiCs molecules separated by a=400 nm; cP/a 3 = 72 kHz. 



the dependence of D and J on the angle 9 for a fixed field magnitude. Note that both D 
and J vanish at 9* = arccos(l/A/3). 

Figure [2b shows that the matrix elements D and J have the same sign, independent of 
9. As discussed in section the sign of the effective mass is always opposite to the sign of 
J, and, consequently, the sign of D. The constant D determines the dynamical interaction 
potential ffTBl) . Due to the linearity of the Schrodinger equation, a positive potential is 
attractive for a particle with negative mass and a negative potential is attractive for a 
particle with positive mass [59J. Therefore, the dynamical interaction in this system is 
always attractive [01]. The kinematic interaction is always repulsive because it arises as a 
consequence of the hard-core boson nature of the excitation operators. 

Figure [2] shows that the values of D and J can be tuned in a wide range of magnitudes. 
For the chosen example of LiCs molecules on an optical lattice with a = 400 nm, D can 
be varied in the interval from -80 kHz to 40 kHz and J from -50 kHz to 25 kHz. Their 
relative magnitudes determine the relative contributions of the kinematic and dynamical 
interactions. By changing the ratio between D and J, one can explore the interaction 
regimes dominated by the different interactions. In particular, prevailing attraction may 
result in the formation of bound two-exciton complexes known as biexcitons [01]. They 
form in the low dimensional systems when \D\ > 2\J\ [65J. In ID they appear as a single 
state band split form the continuum of the two-particle states. As a consequence of the 
correlation between the signs of D and the effective mass, the biexciton state appears below 
(above) the two-particle continuum for negative (positive) D, and can be continuously 
tuned between these two positions by varying the angle between the electric field and the 
intermolecular axis. The properties of biexcitons can be engineered by varying the ratio 
D/2J by tuning the magnitude of the electric field. In particular, for an ensemble of LiCs 
molecules, the biexciton begins to appear at Sf ~ 3.6 kV/cm, when D = 2 J. With the 
increase of the electric field its binding energy (and, correspondingly, its splitting from the 
continuum states) increases, and the wave function shrinks, so that at D ^> 2 J a biexciton 
is a strongly correlated state of two molecular excitations separated by one or two lattice 
constants. The formation of a biexciton thus resembles the association of a molecule by 
combining two atoms. 
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When D = 0, the dynamical interaction vanishes, and the kinematic interaction dom- 
inates. It should be mentioned that the effect of the kinematic interactions on Frenkel 
exciton dynamics has not yet been observed in experiments. To this end, it would be 
useful to find a mechanism for tuning the kinematic interaction as well. We explore this 
in the following section. 

4.3 Suppression of kinematic interaction 

Although the kinematic interaction is inherent to molecular crystals as an intrinsic conse- 
quence of the exciton operator statistics, it may be possible to generate excitations that 
do not experience kinematic interactions. To find the conditions for such excitations, we 
neglect the dynamical interactions and look for solutions of the two-particle Schrodinger 
equation 

H m \*K) = E k \*k) (25) 

in the form 

\^k) = Y,CK(k)P\K/2 + k)P\K/2-k)\v aC ), (26) 
k 

where k = (k± — k 2 )/2 is the relative wave vector of two interacting excitons, and K = 
k\ + k 2 . The expansion coefficients in Eq. ( 12 6 p satisfy the following equation 

(e K {k)-E K y K {k) = (27) 

k 

where ex{k) = E(K/2 — k) + E(K/2 + k) is the total energy of two interacting excitons. The 
two-particle amplitude in the site representation, which depends on the relative distance 
r = n — m between two molecular excitations, is Ck{t) = e lkr CK{k). Since Cx(r = 

k 

0) = 0, the right-hand side of Eq. (1271) vanishes. Therefore, Eq. (1271) can be satisfied only 
if 

e K (k) = const = E K . (28) 
The corresponding wave function is (A^ — l)-time degenerate (n-degeneracy): 

l*i?> = 4ffT. B K{n)e lkn P\K/2 + k)P\K/2-k)\^ C ), (29) 
VM k 

where the quantum number n determines the fixed distance r = na between the excited 
molecules: n — 1,2, Af — 1. These states describe two correlated excitations that do not 
experience the kinematic interaction. 

The key requirement ( 12"S|) - the k- independence of the two-particle energy ex{k) - 
is however not easy to satisfy. One possibility is to consider a specific choice of K. In 
the nearest neighbor approximation, when E(q) = Eq + 2 J cos aq, the two-particle energy 
written in terms of K and k is Sxik) = 2£ , +4Jcos(aA;) cos(aA^/2). It reduces to a constant 
if K = it I a. Thus, with the accuracy up to corrections coming from the interactions 
beyond the nearest neighbors and quadrupole interactions, such a pair immune to the 
kinematic interaction can be produced starting from a pair of excitons with qi + q 2 = tt/cl. 
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Combined with the possibility to control the dynamical interaction strength by the electric 
field magnitude, this provides a system in which both attractive and repulsive interactions 
between quasiparticles can be independently tuned. 

Another possibility to produce a two-particle state with the energy independent of k 
can be realized using the three- level structure of the first rotationally excited state of x £ 
molecules (see Figured^). At low electric fields, when the energy separation between the 
levels \N — 1, M — 0) and \N = 1, M = ±1) is small, the different excitation branches 
may be mixed, which corresponds to configuration mixing in solid state molecular crystals 
|53j . Due to the symmetry properties of the molecular system considered here, the state 
with M = is decoupled from the states with M = ±1 if the electric field is parallel or 
perpendicular to the molecular ensemble. For other angles 8, the collective excitations give 
rise to three exciton modes which we denote a, (3 and 7, containing contributions from all 
three molecular transitions [59] : 

^W*) = £ (30) 

Af =-1,0,1 

where Pl{k) = {l/y/77) £ n e lkn Pl. m are the exciton operators corresponding to the molec- 
ular transition to the excited state with projection M = M^. Operators (13"U|) diagonalize 
the four-level excitonic Hamiltonian yielding 

H^ } = E E p (k)P}(h)P p (k). (31) 

fc;p=o,/3,7 

At high electric field, the 7-mode is split from the other two modes, leading to the iso- 
lated one-band exciton discussed throughout this article. However, the other two branches 
correspond to the excitations of the degenerate molecular states and remain mixed. Their 
k = state can be accessed by microwave field with circular polarization. The non-zero 
k states can be probed by Raman transitions combining photons with linear and circular 
polarization. In the limit of high electric fields, the energies of these excitons can be written 
as 

E a (k; 6) = E x {£ f ) + 2\a N=Q (£ f )\ 2 j(k) (cos 2 9 - 2/3) , 

V 7 (32) 

Ep{k) = E 1 (S f ) + 2\a N=0 (£ f )\ 2 j(k)/3, 

where E\(£f) is the transition energy between the field-dressed ground state and the field- 
dressed excited state with M = ±1; j(k) = (d 2 /a 3 ) Y^[=i cos(akm)/m 3 , and ajs;{£f) are 
defined in Eq. (fT2"|) . The kinematic interaction in the presence of several branches is more 
complicated, but it can be shown that the condition (128|) also applies to the multi-branch 
problems. It can be seen that, at a particular angle 6 = arccos ^/l/3, 

E a (k; 9) + E p {k) = 2E 1 (€ f ) = const, (33) 

so the kinematic interaction must be absent. Eq. (13"31 is satisfied for arbitrary values of 
K. We note that for the 7-branch, the matrix elements of the dipole-dipole interaction 
vanish at 9 = 6* so the 7-state becomes dispersion-less. However, the matrix elements 
giving rise to the other two exciton states remain non-zero at this angle, see Eqs. (|32|) . It 
is the cancellation of the dispersions of E a {k) and Ep(k) that leads to the suppression of 
the kinematic interaction. 



13 



5 Rotational excitons in a tunable disordered poten- 
tial 

Dynamics of quantum particles in disordered potential has been extensively studied in 
relation to Anderson localization [SH], propagation of light [HH] and particles [70] through 
disordered media, phase transitions between insulating and conducting states [71] , to name 
a few examples. As a result of these studies, many properties of disordered systems, such 
as the role of the density of states or the correlation functions are well understood. At the 
same time, there are still many problems that remain open. A few examples include the 
role of disorder in the transition to the glassy state [72], transition from a superconductor 
to insulator with increasing disorder [73] , and the counterintuitive behavior of conductivity 
in some disordered quasi-crystals [73]. 

In this section we show that ultracold molecules trapped on an optical lattice offer the 
possibility to study quasiparticles in the presence of dynamically tunable disordered po- 
tential. We do not consider the effects of unavoidable natural disorder, such as unoccupied 
lattice sites, lattice potential inhomogeneity or fluctuations of the electric field. Always 
present in experiments with atoms and molecules on optical lattices, the natural disorder 
induces localization and decoherence of excitons leading to homogeneous broadening of 
exciton dispersion curves. We assume that the effects of natural disorder can be reduced 
to a small fraction of the exciton bandwidth. Here, we consider the possibility of applying 
an external disordered potential that could be varied to allow the observation of real-time 
dynamics of disorder-induced phenomena. 

5.1 Optical lattice with tunable disorder 

There are several different models of disorder [75], and some of them can potentially be 
realized in an optical lattice with ultracold molecules. First, suppose that a small fraction 
of molecules trapped on the lattice (host molecules) is replaced with molecules of different 
kind (impurities). This results in substitutional disorder. We assume that all impurity 
molecules are identical. The impurities break the translational symmetry of the system 
and therefore scatter excitons. This scattering is elastic and modifies the direction of 
exciton wave propagation, but not the absolute value of the exciton wave vectors. 

An impurity introduced into the molecular crystal in general differs from the host 
molecules by the molecular transition energy, Ei, and by the dipole moment, di, which 
modifies the dipole-dipole coupling strength. In the presence of a single impurity at the 
lattice site n = 0, we can write the total Hamiltonian of the system as 

#tot = H cxc + VoP^Po + £ A J(n) (PlPo + P^P n ) (34) 

where Vq = Ei — Eq and A J(r) is the difference between the host-host and host-impurity 
excitation transfer constants. The exciton-impurity interaction can thus be described as a 
sum of a delta-function potential with strength V and a perturbation due to the difference 
in the dipole moments of the host and impurity molecules. If the matrix of the operator 
( 134)) is evaluated in the basis of exciton states in the site representation, one finds that Vq 
perturbs the diagonal matrix elements and AJ the off-diagonal matrix elements. Thus, 
the constants Vq and AJ(r) give rise to diagonal and off-diagonal disorder. 
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For a properly chosen mixture of diatomic molecules, the magnitude and the sign of Vq 
can be tuned by an external electric field. For example, Figure [3] shows that this can be 
achieved in an array of LiCs molecules doped with LiRb molecules. Generally, it should be 
possible to tune Vq from a negative value to a positive value in a mixture of 1 S diatomic 
molecules AB and CD, when the dipole moment of AB is greater and the rotational 
constant of AB is smaller. 




E f (kV/cm) 



Figure 3: Excitation energies of non-interacting molecules in an electric field for transitions 
\N = 0, M = 0) \N = 1, M) with M = (upper curve) and M = ±1 (lower curve) vs 
electric field for LiCs and LiRb. 

It may also be possible to realize a system with tunable diagonal disorder in a single- 
species ensemble of molecules by applying intense laser beams focused on a small part 
of the lattice. Molecules at the focus of the beam must experience larger ac Stark shifts 
than the molecules outside the beam focus. These molecules are equivalent to impurities 
because their energy level structure is different. The degree of detuning Vq can be adjusted 
by varying the laser field strength. For small fields, it may be possible to neglect the off- 
diagonal disorder terms, which suggests a unique possibility to differentiate between the 
effects of diagonal and off-diagonal disorder. 

5.2 Resonant enhancement of exciton— impurity scattering 

The properties of waves with arbitrary dispersion in the presence of local defects have been 
studied in Ref. [76]. Similar considerations were later applied to excitons in a crystal with 
substitutional disorder in Ref. [77] , where the exciton-impurity scattering cross section, 
which characterizes their interaction, was derived. The exciton-impurity interactions lead 
to the appearance of bound exciton states, which can capture excitons in solid state crystals 
|78j . The constants Vq and A J in Eq. f )3^|) determine the character of the exciton-impurity 
interactions in natural solids. 

We first consider an ensemble of molecules driven by several focused laser beams with 
the same strength, so that A J = 0, and all impurities are characterized by the same value 
of Vq. This system may allow for the possibility to explore the dependence of exciton- 
impurity scattering cross sections not only on the exciton wave vector k, but also on the 
impurity scattering strength Vq, which is not possible in conventional solids. In particular, 
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as we show below, the scattering cross section can be resonantly enhanced, if the potential 
produces a shallow bound state at small values of Vq [59J. 

For a particle with the parabolic dispersion E = h 2 q 2 /2m* in a d- dimensional delta- 
like potential tig = a d V 5(r) the bound state is split from the continuum states by the 
energy Eb, which depends on the dimensionality of the system. The binding energy can be 
obtained from the following equations: 
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where E\ oc = h 2 ir 2 /2\m*\a 2 is the localization energy of a particle with the mass \m*\ in a 
region of dimension a. Note that these equations have solutions only for negative (positive) 
Vo and positive (negative) mass m*. 

Figure UK shows the behavior of the local states as functions of the dimensionless po- 
tential strength |Vo|/-Ei oc - In ID and 2D, an attractive delta-like potential always produces 
a bound state, and in 3D only starting from a finite value: |Vo| > 2E\ oc j-n. Resonant 
scattering may play an important role if the resonant enhancement of the scattering cross 
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section at Eb — > is reached at a finite value of Vq. For instance, in ID E t 
at vanishing Vo, and the shallowing of the bound state is accompanied by the vanishing of 
the scattering potential itself. In turn, in 3D one can expect a resonant enhancement of 
the scattering cross section for the potential \Vq\ ~ 2E\ oc /ix = ixTi /|m*|a 2 . Near this po- 
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tential strength, E t 

produces a shallow bound state, whose energy tends exponentially to zero with vanishing 
impurity potential strength: E^^ ps E\ oc exp[4Eioc/7rVosgn(m*)] (we recall that Vo and 
m* should have different signs for a bound state to appear). This exponential dependence 
leads to efficient resonant scattering for all < Vo < Q.2E\ oc : as the potential strength in 
2D exponentially exceeds the kinetic energy of the scattered wave, then, in contrast to the 
ID case, the resonant scattering cannot be considered as weak, even at vanishingly small 
Vq. At the same time, this means that for | Vo | 0.2_Ei oc the local state is so close to 
zero, that the potential scatters resonantly only excitons with k = 0, which in fact do not 
propagate. 

For small wave vector excitons, E(k) = (Eo + 2J) — Ja 2 k 2 , so that fx 2 jm* = —2Ja 2 , and 
-Eioc = tt 2 J. Borrowing the results from Ref. [79], we write the scattering cross sections in 
2D and 3D for ak 1 as explicit functions of the potential strength Vq, and express them as 
functions of the bound state energy Eb and the kinetic energy T{k) = h 2 k 2 /2|m*| = \ J\(ak) 2 
(these expressions correspond to the limit when both T(k) and Eb are much smaller than 
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Figure 4: (a) Dependence of the bound state energy on potential strength V in 
different dimensions. (b,c) The exciton-impurity scattering cross sections for 2D (b) and 
3D (c) as functions of the potential strength for different values of ak (shown near each 
line), (d) Different propagation regimes in 2D (top) and 3D (bottom) for 1% of impurities 
(thick solid line) and 0.1% of impurities (thin dashed line). 
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Figures Hb andHfc show the scattering cross sections for 2D and 3D (note the divergence 
of cr 3D at k — y 0, Vo — > 2Ei oc /n). These results demonstrate that by tuning the potential 
strength from zero to the critical value (V^ < 0.1E loc and = 2E loc /n), one can 

vary the value of the scattering cross section for small wave vector excitons by many orders 
of magnitude. 

The single-impurity scattering cross sections determine the regimes of exciton propa- 
gation in the lattice with a small admixture of impurities. According to the Ioffe-Regel 
criterion [80J, excitons with the wavelength A = 2n/k are strongly localized when \/2n > I 
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(I is the elastic mean free path). When X/2tc <C I excitons have plane-wave- like charac- 
ter, though, in some cases, weak localization of excitons can be achieved [51] . Excitons 
propagate without scattering when I is on the order of the lattice size. 

For the ensemble of polar molecules trapped on an optical lattice, the elastic mean free 
path of excitons / ~ 1/cxc for a given concentration c of impurities can be dynamically 
tuned by varying the strength of the electric field that modifies the disordered potential 
and, consequently, a. This may allow a possibility to transfer excitons dynamically from 
the regime of ballistic propagation to the regimes of weak or strong localization. Figure HH 
shows the (Vo, /c)-diagram of different propagation regimes for 2D (top) and 3D (bottom) 
lattices with the concentration of impurities 1%. The solid lines correspond to the condition 
I = A/2-7T. Below these lines, the exciton wave length A exceeds 2irl, which corresponds to 
the Anderson localization regime. Above the lines, the states are in general delocalized. 
The dashed line is plotted for the impurity concentration 0.1%. The diagrams are not 
complete without another important parameter, the phase-breaking length l^. This length 
accounts for inelastic scattering processes. The interference effects leading to weak and 
strong localization are only possible if / <C 1$ [81J. The phase-breaking length should be 
calculated for a given realization of the experimental system with the account of the major 
loss channels for rotational excitons. 

5.3 Disorder correlations and localization— derealization crossover 

Consider now a lattice with multiple impurities formed by molecules of a different type. 
This leads to significant values of A J ^ 0. Quantum particles in the presence of a random 
distribution of scattering centers undergo coherent localization, and in ID all states are 
exponentially localized even in presence of disorder of arbitrarily small magnitude |82j . 
However, when the disorder potential exhibits short-range correlations, particular states 
may become delocalized [S3]. Delocalized states may even form a continuous band, if the 
correlations are long-range, so that a mobility edge exists between localized and delocalized 
states [M]. Studies of correlation-induced derealization of quantum states are important 
for understanding quantum transport in disordered systems. 

The substitutional disorder introduces both the diagonal (Vo) and off-diagonal (AJ) 
disorder. As a consequence, Vq and A J are correlated, and the spectrum of a ID disordered 
lattice with a two- molecule mixture must always contain one delocalized state [85] . The 
energy of this state is determined by the relation between Vq and AJ. The derealization 
occurs because the diagonal and off-diagonal perturbations compensate one another. 

As demonstrated by Figure [3J Vq can be tuned by shifting the rotational levels of host 
and impurity molecules simultaneously using a static electric field. In particular, at the 
field corresponding to the encircled region, Vq vanishes, while AJ remains finite. Tuning 
the electric field around this value allows for the possibility to exploit both the positive 
and negative values of Vq. Choosing the appropriate value of Vq can be used to induce the 
derealization of any eigenstate of the system. For states near the origin of the Brillouin 
zone, derealization occurs at Vq ~ — 4A J [55] . Figure shows the evolution of a particular 
eigenstate with V . At V = 0, the state is localized due to non-zero AJ and exhibits the 
characteristic exponential profile. At Vq = 22 kHz, the diagonal disorder compensates the 
effects of the off-diagonal disorder, and the state becomes delocalized. Tuning Vq further 
results in the localization of this eigenstate. 
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Figure 5: Probability density |\&(x)| 2 describing an exciton near the top of the energy 
spectrum for a ID array of 1000 LiCs molecules with 10% homogeneously and randomly 
distributed LiRb impurities. Panels correspond to different values of V : (a) V = 0, (b) 
Vo/h = 22 kHz, and (c) Vo/h = 100 kHz. The difference of the dipole moments of LiCs 
and LiRb molecules leads to the value A J = -6.89 kHz. Figure is taken from Ref . |59j . 



6 Quantum simulation of Holstein polaron model 

Rotational Frenkel excitons can be used for quantum simulation of model Hamiltonians 
that cannot be engineered with atoms, molecules or photons as probe particles. Here, 
we discuss an important example illustrating the possibility of engineering the Holstein 
polaron model with collective excitations of molecules on an optical lattice. 

Polaron is a quasiparticle that describes an electron in a crystal lattice dressed by lattice 
phonons. The interaction properties of polarons are currently researched in an effort to 
understand the mechanism of high-T c superconductivity [86J and quantum transport in 
open quantum systems [87]. Quasiparticles similar to polarons can be created by placing 
an impurity in a Fermi degenerate gas of ultracold atoms, as demonstrated in several 
recent experiments [21]. The impurity, produced by changing the internal state of one 
of the ultracold atoms, can be coupled to the Fermi sea via a Feshbach resonance. This 
gives rise to Fermi polarons. However, lattice phonons are bosons. Therefore, a better 
model of electrons in solid state crystals should be based on coupling a probe particle to 
bosons. While this may be achieved by placing an impurity in a Bose-Einstein condensate 
of ultracold atoms, no such experiments have been reported to date. 

In a recent study [63], Herrera and Krems showed that the dipole-dipole interactions 
between polar molecules on an optical lattice can be exploited to engineer controllable 
couplings between excitons and lattice phonons. In this system, the phonons are associated 
with the oscillatory motion of molecules in the lattice potential 

V = J2Mu 2 5r 2 n /2 + J2U g /\n - m| 3 (37) 

n n,m 

where 5r n is a small deviation from the equilibrium position of the molecule in site n, M 
is the mass of the molecules, ujq is the trapping frequency, and U g = R 3 (gg\Vdd(R)\gg) ■ 
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The first term in Eq. (|37[) describes uncoupled oscillations of molecules in their respective 
lattice sites and depends on the intensity of the trapping laser that determines the trapping 
frequency ljq. The second term accounts for the collective motion by coupling molecules in 
different sites. By introducing the operators aj,(fc) and a u (k) for the phonon mode v with 
wave vector k, we can write the phonon Hamiltonian as 



(3f 



where 



hu ph (k) = w , /l + {12U g /a 5 Mu%) ^ (1 - cosmk)/m 5 . 

V m>0 



(39) 



The phonon spectrum is gapped, as co(k — > 0) — > w , and resembles that of optical phonons 
in solid state crystals. 

The exciton-phonon interaction is obtained by expanding the matrix elements of the 
dipole-dipole interaction in a Taylor series to yield: 



#exc- P h= E g u F (k)[ai{k) + a v (k))P^P n + 

k,u;n 



(40) 



where g v F {k) w F/a^Mu , g v j(k) J/ay/Mu , F = {eg\V dd (a)\eg) - (^|Fdd(a)|^), and 
J is the excitation transfer matrix element defined in Eq. (JHJ). The first term in Eq. ([4*0]) 
describes the phonon-modulated transition energies, the second - the phonon-modulated 
excitation transfer. 
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Figure 6: F and J as functions of the electric field magnitude at a fixed angle 9 = 90° (a) 
and of the angle 9 at dSf = B e (b). Calculations are for a ID ensemble of LiCs molecules 
separated by a=400 nm; d 2 /a 3 = 72 kHz. This is a modified version of the figure taken 
from Ref . [63] . 



The exciton-phonon coupling constants gp and gj depend on the molecular states \g) 
and |e), which can be tuned by an external electric field. The dependence of the constants 
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F and J on the electric field magnitude Ej and direction 9 is shown in Figure El In the 
limit gp 3> g,j the Hamiltonian (HUj) describes the Holstein polaron model [88]; this limit 
corresponds to large values of £f. In the opposite limit, when gp <C gj, it corresponds to 
the model of particle-boson coupling by Su, Schrieffer and Heeger [89] ; this is the limit of 
small £f. 




Time (|is) 

Figure 7: Excitation energy transfer in a ID array of five LiCs molecules separated by 400 
nm in a dc electric field perpendicular to the array, (a) Evolution of the excitation proba- 
bility for the first molecule of the array, when no phonons are present. The revivals occur 
due to the finite size of the molecular array and the reflection of the exciton wavepacket 
from the edge of the ensemble, (b) The same as in (a), but with phonons in an optical 
lattice with trapping frequency = uq/2tc varying in time as indicated in the inset. The 
field strength is 0.5 kV/cm. The figure is taken from Ref . |63] 

The exciton-phonon coupling for molecules in an optical lattice can also be controlled 
by varying the trapping frequency Uq, which is proportional to the intensity of the trapping 
laser and enters both of the coupling constants as 1/ ^JuJq. Increasing the exciton-phonon 
coupling strength results in an increase of the effective mass of the polaron. In the limit 
of strong coupling, polarons are very massive and localized. The coherent propagation is 
impeded and polarons can only propagate by random hopping. By decreasing ujq and there- 
fore increasing the exciton-phonon coupling, it may be possible to induce self-localization 
of excitons. The dynamics of an excitation initially produced on a single molecule is illus- 
trated in Figure [7J These results demonstrate that polar molecules trapped on an optical 
lattice can be used as a dynamically controllable simulator of open quantum systems. 

7 Conclusions 

The development of experimental techniques for cooling, trapping and controlling atoms 
and molecules has opened exciting opportunities for new studies of quantum many-body 
systems. Ultracold atoms and molecules are extensively researched as paradigm systems for 
engineering novel states of quantum matter and for quantum simulation of model Hamil- 
tonians used in condensed matter physics. Most of these studies use ultracold atoms or 
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molecules as probe particles. Here, we discuss the possibility of using collective excita- 
tions in an ensemble of molecules trapped on an optical lattice as a probe in a quantum 
simulation experiment. 

We have shown that rotational excitations of polar molecules on an optical lattice may 
lead to the formation of collective many-body excitations - rotational Frenkel excitons. 
These excitons have unique properties that allow tuning the linear and non-linear exciton 
interactions by modifying the rotational structure of ultracold molecules by an external 
electric field. We suggest that this can be exploited for the study of new regimes of Frenkel 
exciton physics and the dynamics of quantum localization in disordered systems. We 
also suggest that rotational Frenkel excitons can be used for quantum simulation of the 
Holstein polaron model. This offers interesting possibilities to study quantum transport in 
open quantum systems with controllable interactions with the environment. In particular, 
the finite size and tunable properties of the phonon bath for molecules on an optical lattice 
suggests the possibility of exploring the transition from a non-Markovian to Markovian 
environment. To observe rotational excitons, one can measure the populations of the 
rotational states at different lattice sites. As described in Ref . |90j , this can be achieved by 
applying a gradient of an electric field and detecting resonant transitions from Stark-shifted 
molecular levels. 
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